#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>
#include <time.h>

// Multiply a csr matrix with a vector x, and get the resulting vector y ,sum
// use kekan sum
void spmv(int n, int* row_ptr, int* col_idx, double* val, double* x, double* y)
{
    for (int i = 0; i < n; i++) {
        y[i]     = 0.0;
        double c = 0.0;
        for (int j = row_ptr[i]; j < row_ptr[i + 1]; j++) {
            double num = val[j] * x[col_idx[j]];
            double z   = num - c;
            double t   = y[i] + z;
            c          = (t - y[i]) - z;
            y[i]       = t;
        }
    }
}

// Calculate the 2-norm of a vector ,sum use kekan sum
double vec2norm(double* x, int n)
{
    double sum = 0.0;
    double c   = 0.0;
    for (int i = 0; i < n; i++) {
        double num = x[i] * x[i];
        double z   = num - c;
        double t   = sum + z;
        c          = (t - sum) - z;
        sum        = t;
    }

    return sqrt(sum);
}

void check_correctness(int n, int* row_ptr, int* col_idx, double* val, double* x, double* b)
{
    double* b_new   = (double*)malloc(sizeof(double) * n);
    double* check_b = (double*)malloc(sizeof(double) * n);

    spmv(n, row_ptr, col_idx, val, x, b_new); // b_new = Ax
    for (int i = 0; i < n; i++) {
        check_b[i] = b_new[i] - b[i]; // Ax-b
    }
    double answer = vec2norm(check_b, n) / vec2norm(b, n);

    fprintf(stdout, "Check || b - Ax || 2 / || b || 2 =  %12.6e\n", answer);

    free(b_new);
    free(check_b);
}

// return now time
double GetCurrentTime()
{
    struct timeval time;
    gettimeofday(&time, NULL);
    return time.tv_sec * 1000.0 + time.tv_usec / 1000.0;
}